Fourier-synthesis approach for static charge-density reconstruction from theoretical structure factors of CaB6

A novel type of Fourier-synthesis approach is reported for determining electron-density distributions and their Laplacians from static structure factors of CaB6. The approach relies on mathematical weighting functions to yield a data set, reproducing all characteristic chemical bonding features of the original quantum-chemically calculated distributions.


Introduction
As a consequence of the Hohenberg-Kohn theorems (Hohenberg & Kohn, 1964), the electron density (ED) of a chemical system is known to constitute its most fundamental observable property. Historically, in early X-ray diffraction experiments, only the locations of the atoms were identified from the dominating local maxima of the ED, which can be successfully modelled in the independent atom model (IAM), i.e. using a superposition of free atom form factors to model the structure factor for each reflection. With the improvement of experimental setups, the reconstruction of ED distributions from diffraction experiments came into reach, and methodological developments were devised for reliable extraction of this information from the observed diffraction intensities of elastic X-ray scattering experiments. Two fundamentally different strategies were followed for this purpose (Waser & Schomaker, 1953). One possibility consists of setting up a position-space model for the ED distribution and fitting it to the experimentally observed intensities. This was realized by replacing the superposition of spherically symmetric ED distributions in the IAM by a multipolar expansion, with the multipoles at each atomic site being consistent with the corresponding site symmetry. The most prominent of these methods used nowadays is the Hansen-Coppens (HC) multipole model (Hansen & Coppens, 1978). On the other hand, it was recognized early on that, in principle, the ED distribution can be reconstructed by 'Fourier inversion' of structure factors (Bragg & West, 1930). If this procedure is terminated before the remaining members of the Fourier series display negligibly small structure-factor amplitudes, series truncation artefacts like spurious local ED maxima and minima can occur, which adversely affect chemical interpretation. In the following, the process of Fourier inversion using a limited number of structure factors eventually multiplied by a weighting function will be called Fourier synthesis. It has long been known that multiplication of the static structure factors by the Debye-Waller factor helps to avoid series termination artefacts (Bragg & West, 1930). However, this turns the static ED distribution into a physically different, temperature-dependent dynamically smeared one. In pioneering experimental studies on Fourier synthesis of ED distributions in prototype ionic, polar covalent and covalent bonding situations in NaCl, MgO and diamond, the introduction of artificially enhanced temperature factors was found necessary in vibrationally hard compounds in order to eliminate series truncation artefacts within experimental resolution (Brill et al., 1939(Brill et al., , 1948. Avoiding artificial temperature factors, a recent study reported on Fourier-synthesized dynamically smeared all-ED (all-electron density) distributions for two amino acids free from series truncation effects (Mondal et al., 2012). They were calculated from static structure factors and atomic displacement parameters obtained by HC multipole model fits of temperature-dependent X-ray diffraction data with resolutions [sin()/] max 1.2 Å -1 . Good convergence of the synthesized dynamical all-ED and all-ED Laplacian distributions was shown to require usage of extrapolated structure factors from the fitted multipole model up to cut-offs at [sin()/] max ! 6 Å -1 . Note that this extrapolation to high resolution was necessary because of a combination of very strict numerical convergence conditions with comparably small physical vibrational smearing working against series termination effects. In order to perform Fourier synthesis using only HC dynamical structure factors within experimental resolution, a suitable mathematical weighting function may be invoked.
Notably, even the practice in standard HC studies on experimental EDs to calculate the static ED distributions from fitted atomic multipoles in real space (avoiding Fourier synthesis) corresponds to an implicit extrapolation of the HC model structure factors to infinite resolution (Coppens & Stevens, 1977). An alternative route more consistent with the experimental conditions would be the Fourier synthesis of the ED using only those static structure factors of the HC model that correspond to the experiment resolution. To the best of our knowledge, this has not yet been done.
In the present methodological study, the focus lies on mathematical weighting functions in order to keep the physical interpretation of the original data set, i.e. Fourier synthesis using static ED structure factors and mathematical weighting functions still yields (smoothed) static ED distributions. The effectiveness of a number of mathematical weighting functions on convergence of Fourier synthesis of static valence (val)-EDs and all-EDs and their Laplacians is investigated. One reason for selection of static ED Fourier synthesis is that there is a theoretical foundation for a QTAIM (quantum theory of atoms in molecules) (Bader, 1990) topological study for static EDs. Moreover, the focus on static ED synthesis represents a more challenging task compared with the dynamic one, because series termination effects and artefacts are more prominent there. By converging the Fourier synthesis of static EDs already at experimental resolution with suitable mathematical weighting functions, the extrapolation of structure factors or usage of artificially enhanced Debye-Waller factors can be avoided. The same technique may also work for Fourier synthesis of dynamic EDs, thus avoiding extrapolation of experimental data sets (see above).   (200) and (110) planes with molecular graph; À i (rank, signature) enumerates and classifies the crystallographically independent ED critical points; ED isoline steps at 0.01 e À bohr À3 , ED Laplacian isoline steps of 0.05 e À bohr À5 .
Finally, the HC model fitting is not free from difficulties and biases (Michael & Koritsanszky, 2017), such that the fitted structure factors {F HC } may be different from the initially given ones {F 0 }. It is conceivable that, in certain cases, significantly different ED or ED Laplacian distributions may be obtained compared with those obtained from Fourier synthesis using the initial structure-factor set {F 0 }.
On the methodological side, static EDs in general are available only from a model. In the present case, a quantumchemical DFT (density functional theory) calculation on a crystalline system was employed, and the ED Fourier coefficients (structure factors) can be obtained up to arbitrary resolution. The static ED and its Laplacian calculated from the model correspond to the reference distributions in the subsequent investigations, and the various Fourier-synthesized distributions were evaluated with respect to these references not only at critical points, but also in the whole unit cell and parts of it using norm deviations.
As an example test case, the non-molecular crystal structure of CaB 6 has been chosen (employed experimental structure data in Table 1). The proper reconstruction of its ED features, which are related to a rather evolved chemical bonding scenario, is considered to be methodologically more challenging than for a molecular crystal of a classical organic molecule like, e.g., glycine or urea. The cubic structure [space group Pm " 3 3m, Fig. 1(a)] of binary hexaborides MB 6 (M = Na, K, Rb, Ca, Sr, Ba, Sc, Y, La) is based on a CsCl-like arrangement of metal atoms and barycentres of sixfold interconnected B 6 octahedra. The octahedral B 6 clusters display six short interoctahedral (exohedral) and 12 longer intra-octahedral (endoendohedral) B-B contacts with substantially different bonding character. Conceptually, the exohedral contacts represent single two-centre (2c) bonds, while the endohedral ones are described by a mixture of 2c and three-centre (3c) bonding character (Longuet-Higgins & Roberts, 1954). A recent theoretical study focuses on the description of chemical bonding for binary hexaboride compounds using positionspace bonding indicators, namely the ED distribution (r), the ED Laplacian r 2 (r), as well as the pair density based electron localizability indicator distribution (ELI-D), and 2c and 3c delocalization indices (DIs) between atomic regions (Bö rrnert et al., 2013) (the first author of the present paper, Dr Carina Bergner, was previously known as Dr Carina Bö rrnert). The resulting QTAIM effective charge Ca 1.52+ is consistent with formal charge assignment Ca 2+ obtained on the basis of ELIBON (ELI-D-based oxidation number). This is in accordance with the saturation of electronic demand of the boron framework composed of 6-connected closo B 6 2À clusters by charge transfer of two electrons from Ca according to the Lipscomb-Wade rules (Lipscomb, 1979;Wade, 1971). Two-and three-centre DI analysis reveals a 2c-2-electron character of the exohedral B-B bonds, and -in agreement with ELI-D topology -a mixed 2c + 3c character of endohedral B-B bonding. Endohedral 2c DIs of about 0.60, and 3c ones of about 0.20 are found for electronically saturated hexaborides of the alkaline-earth metals. These values are related to the total amount of seven endohedral bonds expected according to Lipscomb-Wade rules (Bö rrnert et al., 2013).
The pair density based bonding indicators ELI-D and DIs are still not accessible in purely experimental studies of all types of crystalline compounds. The very useful X-ray constrained wavefunction calculation method, which delivers this information at a semi-experimental level, is nowadays still technically limited to molecular crystals (Genoni & Jayatilaka, 2021). Using diffraction experiments, the ED distribution can be reliably reconstructed from structure factors. Certain specific features of chemical bonding in hexaborides can already be obtained from ED-based quantities, e.g. QTAIM atomic charges and charge-density concentrations [r 2 (r c ) < 0] at bond and ring critical points r c . They fit the generally accepted bonding picture. In the quantum-chemical calculations on CaB 6 , comparably high ED is found within the B 6 skeleton, i.e. edges and faces (Mebs et al., 2011), and lower values in the octahedron centre and interstitial regions [ Fig.  1(b)]. The topological features of the ED distribution, namely the critical points r c obeying r(r c ) = 0, have been determined within the unit cell and characterized by their curvatures (Bö rrnert et al., 2013;Bö rrnert, 2013). The set of all critical points of types (rank, signature) within the unit cell [ Fig. 1(b)] satisfies the Poincaré -Hopf relationship (Zou & Bader, 1994). The distribution of the ED Laplacian [ Fig. 1(c)] contains additional information about particular chemical bonding features in CaB 6 . Values of r 2 < 0 at the bond critical points (b.c.p.) À 4 [d(B-B) endo = 1.751 Å ] and À 5 [d(B-B) exo = 1.676 Å ] between the boron atoms indicate a charge concentration which is characteristic for covalent bonding within the boron framework. Apart from the Ca-atom locations, the inter-octahedral region of the CaB 6 structure is characterized by low values of . A b.c.p. À 3 is found between the octahedral face and the Ca atom [ Fig. 1(b)]. The positive value of r 2 at and close to this b.c.p. represents charge depletion and indicates a predominantly ionic interaction between the metal atoms and the boron framework. Negative values of r 2 not only along the endohedral B-B bond path and b.c.p. À 4 , but also around the B 6 clusters, and especially at the ring critical points (r.c.p.) À 6 slightly above the octahedral triangular face midpoints (Bader & Legare, 1992), mark the whole endohedral boron skeleton (i.e. B-B bonds and B-B-B triangles) as a region of charge-density concentration (r 2 < 0). Chargedensity depletion is found for the interior part of the octahedron, with the highest depletion found at the cage critical point (c.c.p.) À 8 at the octahedron centre. Summarizing, it is evident that for CaB 6 ED characteristics qualitatively feature certain results of pair density based bonding indicators important for chemical understanding of this type of compound. The purpose of this work is to investigate, on the basis of the calculated ED and ED Laplacian distribution for CaB 6 (structure parameters, Table 1) from a periodic DFT calculation, whether its decisive chemical bonding features can be extracted from a data set of static structure factors with a given resolution 1 2 H max = [sin()/] max . Within this analysis, reconstructed distributions of and r 2 for CaB 6, obtained by different Fourier-synthesis methods, are compared with each other and with the original distributions obtained from DFT calculation at the experimental geometry. All structure factors were calculated from the DFT wavefunction (see Appendix A). Comparison of the relevant chemical bonding features of the quantum-chemically calculated static reference ED with those obtained from a suitable Fourier synthesis at experimentally accessible resolution [sin()/] max can give an answer to the question of whether that resolution is sufficient to yield the chemical bonding information content of the original distribution. This means a suitable Fourier-synthesis procedure may be a valuable tool to evaluate the ideally extractable information content of the original ED from a certain measured or measurable data set with limited resolution.

ED and ED Laplacian distributions from theoretical structure factors by the method(exponent) type of Fouriersynthesis approach
The ED (r) of a crystalline system is physically observable. Diffraction measurements on crystalline samples yield diffraction intensities I(H) at reciprocal-lattice positions H ¼ ha Ã þ kb Ã þ lc Ã , which are related to the squared modulus of the product of the corresponding static structure factors F(H) and the Debye-Waller factor. These materialand structure-dependent quantities are obtained by fitting intensities of a coherent elastic scattering experiment, e.g. by means of the HC pseudo-atom model (Hansen & Coppens, 1978;Gillet & Koritsanszky, 2012).
The ED structure factors are defined by a Fourier integral of the corresponding ED distribution (r) in the unit cell (u.c.), the static or the dynamically smeared one. For the purpose of a subsequent QTAIM analysis, we focus on static EDs and structure factors F(H) in the following, Mathematical inversion of this formula leads to the well known Fourier-series representation of the ED (Waser & Schomaker, 1953): Exact reconstruction of the original ED via Fourier backtransformation is impossible in practice due to the necessity of infinite structure-factor sets; all procedures working with such incomplete back-transformations will be called Fourier synthesis in the following, The incompleteness of Fourier series S n (r) often leads to negative ED values, spurious non-nuclear maxima (NNMs) and other errors, e.g. shifted nuclear ED maxima, broadened nuclear ED 'peaks' and ED ripples (Altomare et al., 2008), seen in the raw Fourier-synthesized [equation (3)] EDs [Fig. 2(a)]; all of them are typically classified as series termination errors. In the present Fourier-synthesis approach for ED and ED Laplacian reconstruction for QTAIM analysis, we tend to distinguish between series truncation artefacts (negative ED values and NNMs not present in the reference ED) and so-called (normal) series termination errors. The latter correspond to all other deviation types of a limitedresolution Fourier synthesis from the mathematically exactly [at infinite resolution, equation (2)] back-transformed ED. An enumeration of deviation types is typically dependent on context, such that for the present QTAIM type of ED analysis, additional series termination errors were even more relevant than those mentioned before, namely shifts of positions of critical points, deviations of ED and ED Laplacian values at critical points, and ED Laplacian ripples. The distinction between artefacts and errors reflects the evaluation strategy adopted in the present approach. Fourier-synthesis results with ED features classified as series termination artefacts are considered of minor quality and are usually dropped. In contrast, Fourier-synthesis results displaying (normal) series termination errors (in the present definition) are principally accepted while working systematically to decrease them (see below). The classification as series termination artefact is clearly possible for the negative ED values, which are not allowed for real EDs; the situation is more delicate for NNM occurrence. NNMs are principally allowed features for real ED distributions, though rather seldom occur. Their occurrence during Fourier synthesis of an unknown original ED should be carefully checked in order to either establish their reality or identify them as artefacts. This would be related to potential future applications as a stand-alone method and not for the present study, where they are clearly identified as artefacts, because the reference DFT ED from which the ED structure factors have been obtained does not display NNMs.
Weighting-function methods. S n ðrÞ is mathematically defined by folding the true distribution of (r) with the Dirichlet kernel D n (r). The latter is a cosine series oscillating in sign and value around D n (r) = 0 (Pepinsky, 1952). This is the reason why Fourier synthesis with this kernel ('raw Fourier synthesis') may even yield regions in space where S n (r) < 0. ED distributions obtained from partial sums S n (r) are characterized by Euclidean norm deviations 2 n , which monotonically decay with increasing number n of partial sum elements (L 2 norm) (Carleson, 1966): For a given distribution S n (r), a reduction of the truncation errors can be achieved by modifying the Dirichlet kernel, evoking a gradual down-weighting of partial sum elements, i.e. structure factors, with increasing length jHj ¼ H of the scattering vector. This automatically occurs already by multiplication of the static structure factor with the Debye-Waller factor, which decays as exp(-1 4 B iso H 2 ) (Coppens, 1997). It results from folding the static ED distribution with the thermal motion of the atoms. In order to smooth truncation errors, early works on ED reconstruction either suggested the introduction of an artificial temperature factor (Waser & Schomaker, 1953;Pepinsky, 1952) or extrapolation of the measured data using artificial series members (van Reijen, 1942). Acceleration of the convergence of the series expansion by substituting the strongly oscillating Dirichlet kernel by integral kernels with invariably positive terms offers a purely mathematical method to avoid series truncation artefacts and reduce certain errors without using empirical functions (Pepinsky, 1952). The most prominent of those kernels is the Fejer kernel, which is obtained by Cesaro summation, i.e. calculating the arithmetic average C nmax (r) of partial sums S n (r) defined in equation (3): This leads to a maximum weight of 1.0 for the lowest structure factor F(000) = F(H j=0 ), and a monotonic decrease of the weights C (H j ) of the remaining nmax structure factors F(H j ) according to 1 À n/(nmax+1) along the sequence n = 0, 1, . . . nmax. For this reason, the non-negative Fejer kernel acts as a low-pass filter on the structure factors, which leads to uniform convergence of a Fourier series. As can be seen from the 1D index n, the Fejer method historically is a 1D method, where the Cesaro summation of partial sums of Fourier coefficients finally leads to a simple multiplicative factor C (H j ) for each Fourier coefficient F(H j ) in the Fourier-synthesis formula [equation (6)]: For a given resolution [sin()/] max = 1 2 H max and a corresponding number n(H max )+1 of Fourier coefficients to be included, the size for each C (H j ) factor varies in the range ]0, 1] and is related to the position of the Fourier coefficient in the sequence of coefficients j = 0. . . n(H max ), according to C ðH j Þ ¼ 1 À nðH j Þ=½nðH max Þ þ 1, where j = n(H max ) denotes the final coefficient to be included with non-zero weight. Application of this method to structure-factor sets {F(H j )} organized by three indices h j , k j , l j is a non-unique extension. Such a 3D variant for ED reconstruction has been proposed and applied by Pepinsky (1952). It corresponds to C (H j ) factors composed of a product of the three 1D terms related to h j , k j and l j : All symmetry-related structure factors [H j ] sym belonging to the same reflection class obtain the same weighting factor, which is a common feature of all weighting schemes. As a special feature of all 3D weighting schemes presented herein, symmetry-independent structure factors with the same 1 2 |H j | = sin()/ obtain different weighting factors. Note that this weighting method denoted Pep3D in the following has some intrinsic disadvantage. It selects [equation (7)] a box-shaped region from the available structure factors F hkl , which introduces a certain directional anisotropy. Nevertheless, numerical evaluation of this method has been selected herein also for historical reasons.
Utilizing the degrees of freedom in the weighting schemes denoted as 1D in the following, an ordered string (1D) of structure factors {F(H j )} is formed using the natural sequence of increasing H j = |H j |. So, the ordered string looks like {[H j ] |H| } ord . Not only will all symmetry-related structure factors be given the same weighting factor, but so will all reflections with the same H j in general. This is a necessary consequence of having a uniquely defined sequence of structure factors with increasing H j leading to monotonically decaying sigma factors with H j . In the Fejer spirit, we may now assign increasing position numbers n j to each class [H j ] |H| . F(0, 0, 0) has to be assigned position number n 0 = 0, as it must always obtain a weighting factor of 1. Three different ways to assign position numbers for the remaining list of structure-factor classes [H j ] |H| were investigated.
In the conceptually simplest case (at first sight), we may assign consecutive integer position number values n j = j = 1, 2, . . . , N max to the sequence of N max classes [H j ] |H| . As an example, the first three classes [0,0,0] |H| , [1,0,0] |H| and [1,1,0] |H| in the ordered list obtain position numbers n j = 0, 1 and 2, respectively. This weighting method is denoted Fej_cnt in the following: In the variant denoted Fej_pcl, we take into account the number of structure factors inside each [H j ] |H| class. The procedure uses the consecutive integer counting numbers j of each single structure factor in the complete ordered list. Each of the N max classes [H j ] |H| obtains the position value of the arithmetic mean of the first and the last member's position value, n j = (j first + j last )/2. As an example, the first classes in the list [0,0,0] |H| , [1,0,0] |H| and [1,1,0] |H| obtain position numbers 0, 3.5 and 12.5, respectively. This procedure may be justified considering that an infinitesimal, e.g. orthorhombic, distortion would keep the number of structure factors, but split the [H j ] |H| classes and lead to additional position numbers compared with Fej_cnt. These are consistently taken into account in the Fej_pcl weighting method: Having accepted position numbers for the [H j ] |H| classes being fractional and not necessarily increasing by 1, the introduction of Fej_stl represents only a small additional step. In this weighting scheme there is an infinite number of position number holes between the present ones, keeping the space for additional structure factors (reflections), e.g. in case of symmetry reductions that lead to vanishing of non-primitive lattice translations. This is effectuated by using the H j values themselves as 'position numbers'. Moreover, this 1D scheme is a directly related alternative to the 3D Pepinsky one. Instead of using the product of three 1D terms, the h, k, l values are additively combined according to H j ¼ 1=aðh 2 j þ k 2 j þ l 2 j Þ 1=2 for the cubic structure [H j = 1/d hkl , i.e. the reciprocal interplane distance of lattice planes (h j k j l j )], with 'a' being the lattice parameter: Since 1 2 H j = sin()/, weighting factor Fej1D will linearly decrease with sin()/. Note that all 1D schemes presented herein automatically avoid the problem of selecting structures F hkl in a box-shaped region like in the original Pep3D [equation (7)]. The region defined by H max is always spherical in reciprocal space.
Based on Lanczos's solution of smoothing the first derivative of a 1D Fourier expansion (Lanczos, 1956), Shchedrin & Simonov (1969) suggested an algorithm for the 3D-expanded ED Fourier-synthesis problem. Although the authors were not interested in the ED gradient itself, it was recognized that smoothing the ED gradient simultaneously smooths the contours of the ED itself. Therefore, the Lanczos factors can serve the same purpose as the Fejer factors, the smoothing of the ED to counteract series termination artefacts: This weighting method, herein denoted Lan3D, has been implemented and employed in the following. Note that the special condition [equation (12)] was especially designed to avoid box-type selection of structure factors as occurs in the original Pep3D [equation (7)]. Interestingly, a similar weighting function has also been used for Fourier back-transformation of magnetic structure factors obtained from polarized neutron scattering experiments to calculate experimental magnetization density maps (Shull & Mook, 1966).
In analogy to method Fej_stl, a 1D Lanczos method has been set up, called Lan1D in the following. It is possibly identical to the method mentioned by Strel'tsov et al. (1985), but without the explicit formulas given there. Similar to the 1D Fejer methods, in the Lan1D weighting scheme all structure factors with the same H j obtain the same numerical weightingfactor value, while this is not so in the 3D schemes Pep3D and Lan3D: With these six weighting-function methods for compensation of series termination effects in the Fourier synthesis of the ED distribution at hand, the question remains as to how to achieve a smoothed ED Laplacian distribution?
Weighting-function exponent. While Fourier series in general are convergent, this is not necessarily valid for their derivatives. Caused by the additional factor 4 2 H j 2 , which increases with increasing size (resolution) of the Fourier series faster than the static structure factors decay, the second derivatives of these Fourier series are divergent, As already mentioned, the Lanczos scheme was initially set up to suppress oscillation artefacts in the first derivative of a truncated Fourier series. Clearly, these artefacts are already contained in the synthesized ED distribution itself, and a numerical differentiation of this distribution gives the same artefacts as the analytical one. Therefore, smoothing out series termination artefacts of the ED will simultaneously improve the ED gradient and Laplacian representation as well. The interconnection between the ED distribution, its gradients and its Laplacian distribution is directly employed in the QTAIM topological analysis of the ED distribution. The critical points are located where the ED gradient is zero, and they are characterized by the ED curvatures at this point. Therefore, the ED, the ED gradient and the ED Laplacian have to be synthesized using the same sigma factor p (H j ). This simple conclusion also implies a further degree of freedom in the sigma factor's usage. Since the initial sigma factor derived by Lanczos was for the first derivative, the Laplacian should have the squared sigma factor, both of which cannot be fulfilled simultaneously in the topological analysis of the ED. This not only means that one can use either exponent, 1 or 2, but that one can use any real number p ! 1.0 (and eventually even p > 0). Exploiting this degree of freedom has been found vital for the whole study. To the best of our knowledge, this has never been considered or formulated before in the framework of Fourier synthesis of EDs and ED Laplacians. As a strategy, for each Fourier-synthesis method presented herein, the lowest exponent p has been used, which yields a series of partial sums S n (r) without artefact NNMs at 1 2 H max > 1.1 Å À1 . For the present study shown below, this was achieved with 1.0 p 3.0, with the actual value used depending on the method. The specification of the weighting function used in each case will be abbreviated in the form method(p value), e.g. Lan1D(1.75). Summarizing, for Fourier synthesis the following equations were always employed: In the following, the p superscript of the synthesized distributions S p n ðrÞ [equation (15)] and r 2 S p n ðrÞ [equation (16)] indicating the present method(exponent) (ME) type of approach is omitted for brevity. The effect of the p factors with respect to the convergence of the ED and its Laplacian distribution with increasing resolution 1 2 H max is investigated using the example of CaB 6 .
Although the application of Fourier synthesis for ED reconstruction and even for its Laplacian has been reviewed (Tsirelson & Ozerov, 1996), so far, the approach itself (inde-pendent from experimental uncertainties) has never been systematically investigated with respect to convergence behaviour of QTAIM-related properties. Typically, only the final results obtained from the experimental [sin()/] max are reported, but the behaviour of the final values upon increasing resolution from lower resolutions up to this final one would be a further criterion to judge reliability of the final results. Moreover, as explained above, in a QTAIM study the Fourier syntheses of the ED and the ED Laplacian are conceptually restricted [equations (15), (16)] to employ the same smoothing factors p (H j ). To the best of our knowledge, it has never been explicitly noted in the literature that this has actually been done.
Valence-electron and all-electron Fourier synthesis. A partitioning of all-electron density (all-ED) distributions into additive core-electron (core-ED) and valence-electron (val-ED) distributions has been found useful in chemistry, physics and crystallography, i.e. for the reference DFT all-ED and all-ED Laplacian, one can formulate DFT ðrÞ ¼ core ðrÞ þ val ðrÞ ð 17Þ and r 2 DFT ðrÞ ¼ r 2 core ðrÞ þ r 2 val ðrÞ: The core states of a compound, e.g. B(1s 2 ), Ca(1s 2 , 2s 2 , 2p 6 , 3s 2 , 3p 6 ) for CaB 6 , are assumed to be chemically inert, such that the core-ED can be reasonably approximated by the core-ED from free atoms, and the chemical focus typically lies on characteristic features of the val-ED specific for the respective compound. In the present study, Fourier syntheses of val-and all-EDs and their Laplacian distributions are investigated. They are computed from application of equations (15), (16) using val-ED and all-ED structure factors obtained for the reference DFT wavefunction from the program Elk (see Appendix A). With the subsequent QTAIM analysis in mind, in the present study all-ED and all-ED Laplacians are always employed and evaluated with respect to the corresponding reference DFT all-electron distributions. All-ED and all-ED Laplacian distributions synthesized from all-ED structure factors are denoted S n,tot (r) and r 2 S n,tot (r), respectively, in the following. For synthesized val-ED distributions S n,val (r) and r 2 S n,val (r) only, an additional step is required to construct a special kind of all-ED distribution S n;val ðrÞ and r 2 S n;val ðrÞ by adding the corresponding reference DFT core-ED distributions core (r) and r 2 core (r): S n;val ðrÞ ¼ S n;val ðrÞ þ core ðrÞ ð 19Þ r 2 S n;val ðrÞ ¼ r 2 S n;val ðrÞ þ r 2 core ðrÞ: Calculated deviations (see the next section) of this kind of all-ED [equation (19)] and all-ED Laplacian [equation (20)] with respect to the corresponding reference DFT all-ED [equation (17)] and all-ED Laplacian [equation (18)] can be easily seen, S n;val ðrÞ À DFT ðrÞ ¼ S n;val ðrÞ À val ðrÞ; ð21Þ r 2 S n;val ðrÞ À r 2 DFT ðrÞ ¼ r 2 S n;val ðrÞ À r 2 val ðrÞ; ð22Þ to solely measure the val-ED and val-ED Laplacian deviations with respect to the corresponding reference DFT val-ED and val-ED Laplacian distributions.

Statistical and topological evaluation of Fourier-synthesized EDs and ED Laplacians with respect to deviations from reference DFT-based distributions
For evaluation of obtained ED and ED Laplacian distributions S n (r) and r 2 S n (r) determined with increasing resolutions 1 2 H max by various Fourier-synthesis methods(exponents), several statistical quality and convergence measures have been calculated. They are all based on the norms of the deviations of the obtained distributions from the reference distributions DFT (r) and r 2 DFT (r).
With the ED (r) being a square-integrable function, the sequence of its Fourier-synthesized distributions S n (r) is expected to converge with the number n of structure factors in the Euclidean norm L 2 (Riesz-Fischer theorem, more precisely L m , with 2 m < 1), i.e. L 2 n r u:c: ; S n ð Þ¼ R u:c: A different way to perform Fourier synthesis proceeds via Fejer summation, which can be written as a weighted Fourier sum (see above). This kind of summation is expected to show uniform convergence (Fejer's theorem), and converges already in the L 1 norm, L 1 n r u:c: ; S n ð Þ¼ R u:c: S n ðrÞ À ðrÞ dr ' V u:c: P N vox;u:c: j¼1 S n ðjÞ À ðjÞ ! 0 as n ! 1: For the criterion of 'uniform convergence' to be fulfilled, the maximum norm L 1 has to converge according to S n ðjÞ À ðjÞ ! 0 as n ! 1: ð25Þ Note that uniform convergence always implies Euclidean and point-wise convergence, but not vice versa (Tveito & Winther, 2005). The maximum norm Á max (S n ) describes the maximal deviation between the original and the reconstructed ED in the unit cell. It is a very transparent measure of convergence, although, when calculated on a finite grid, it may be more prone to missing positions with higher deviations than the mean square convergence measure (seen in one case below). Note that the position r of maximal deviation may change with increasing resolution.
The actual computations were performed by discrete summations over an orthogonal equidistant real-space mesh of N vox,u.c. = 158 Â 158 Â 158 grid voxels j (mesh size = 0.02628 Å ) covering the whole CaB 6 unit cell. The all-electron [F tot (000) = 50 e À ] and valence-electron [F val (000) = 20 e À ] density structure factors were calculated from the DFT wavefunction of the full unit cell. The ED and ED Laplacian distributions obtained from Fourier synthesis were evaluated on the uniform grid in the full unit cell using the L n 1 , L n 2 and L n 1 norm deviations with respect to the reference DFT distributions DFT (j) and r 2 DFT (j).
The L n 1 and L n 2 norm deviations per voxel of the reconstructed val-EDs S n,val (r) and their Laplacians r 2 S n,val (r) (replace 'S n;val ' by 'r 2 S n;val '), and all-EDs S n,tot (r) and their Laplacians r 2 S n;tot ðrÞ (replace 'S n,tot ' by 'r 2 S n;tot ') with respect to the reference DFT all-ED DFT (r) and its Laplacian have been calculated on the unit-cell grid positions r u.c. : L 1 n r u:c: ; S n;tot L 1 n r u:c: ; S n;val L 2 n r u:c: ; S n;tot L 2 n r u:c: ; S n;val In the present study, the focus lies mainly on the critical point reconstruction between the atoms. Therefore, the statistical measures introduced have been computed additionally for the valence region between the atomic core regions. Summation over the voxels of the valence region (N vox,val = 3327147) was performed, skipping the voxels inside the spherical regions around the atomic positions with radii of 2.495 and 0.75 bohr around Ca and B atoms, respectively, which are identified as core regions from ELF/ELI-D (ELF = electron localization function) atomic shell structure studies of free atoms (Kohout & Savin, 1996;Baranov, 2014): L 1 n r val ; S n;val L 2 n r val ; S n;tot research papers L 2 n r val ; S n;val While the theorems on the convergence are valid for the full region of the unit cell (denoted r u.c. ), it is not clear how the statistical deviation measures behave in the valence region (denoted r val ). For all distributions calculated herein L n 1 and L n 2 norm deviations of the same distribution were found to obey the known relations L n 1 > L n 2 and L n 1 < N vox 1/2 L n 2 , as mathematically expected. Moreover, for all distributions calculated herein, it was also found that for each Fouriersynthesis function, the absolute deviations in the valence regions r val were smaller than the corresponding ones in the full regions, L n 1 (r val ) < L n 1 (r u.c. ) and L n 2 (r val ) < L n 2 (r u.c. ), which is caused by the rather different maximal ED values in both regions. In the following, it was considered more convenient to compare relative norm deviations with one another, taking into account the different average values in the respective regions. The relative norm deviations were obtained from division of the norm deviations by the DFT or val 1-and 2-norms of the corresponding reference DFT all-ED and val-ED, respectively, i.e. with m = 1, 2, L m n r u:c: ; S n;tot L m n r val ; S n;tot L m n r u:c: ; S n;val L m n r val ; S n;val where L m r u:c: L m r u:c: ; val ð Þ¼ 1 N vox;u:c: In an analogous way, relative norm deviations for the ED Laplacian distributions were calculated as well [replace in equations (34)-(37) 'S n,DFT ' and 'S n,val ' by 'r 2 S n,DFT ' and 'r 2 S n,val ', respectively, and in equations (34)-(41) ' DFT ' and ' val ' by 'r 2 DFT ' and 'r 2 val ', respectively]. The 1-norm per voxel of the all-ED DFT has a simple meaning: it is the average ED. Exact summation of the CaB 6 all-ED DFT in the whole unit cell [equation (38)] would yield L 1 (r u.c. , DFT ) = 50 e À /V u.c. , i.e. the average all-ED. Exact summation of the all-ED in the ELI-D valence region [equation (39)] would approximately yield L 1 (r val , DFT ) '20 e À /V u.c. , which is mainly connected with ELI-D shell structure representation reproducing not quantitatively exactly the integer atomic shell occupations given by the periodic table of the elements (Kohout & Savin, 1996;Baranov, 2014). Note, while the exact sum (integral) of the ED Laplacian over the unit cell has to be zero due to exact cancellation of positive and negative values, the 1-norm of the ED Laplacian being the sum (integral) of its absolute values is mathematically unbounded from above. The values of the employed L 1 and L 2 norms of the reference DFT ED and its Laplacian are given in Table 2. They all correspond to those obtained by summations on the equidistant grid specified above. Thus, some deviation with respect to the theoretically expected average ED values can be observed, which plays no role in comparison of relative norm deviations of Fouriersynthesized distributions normalized with the same reference DFT ED norm.
Of interest is also the largest absolute deviation of each S n (r) and r 2 S n (r) distribution in the whole unit cell and in the valence region, denoted according to Á max (r u.c. , S n,tot ), Á max (r val , S n,tot ), and Á max (r u.c. , S n,val ), Á max (r val , S n,val ), respectively: Á max r u:c: ; S n;tot À Á ¼ L 1 r u:c: ; S n;tot À Á ¼ sup j2u:c: S n;tot ðjÞ À DFT ðjÞ ð42Þ Á max r u:c: ; S n;val À Á ¼ L 1 r u:c: ; S n;val À Á ¼ sup j2u:c: S n;val ðjÞ À DFT ðjÞ : Thus, the following 12 deviation measures for the synthesized ED distributions have been calculated in the present study: L n m (r u.c. , S n,tot ) rel.t , L n m (r u.c. , S n,val ) rel.v , L n m (r val , S n,tot ) rel.t , L n m (r val , S n,val ) rel.v , with m = 1, 2, and Á max (r u.c. , S n,tot ), Á max (r val , S n,tot ), Á max (r u.c. , S n,val ) and Á max (r val , S n,val  L 1 and L 2 norms (in atomic units) of the reference DFT ED and ED Laplacian distributions ( val , DFT and r 2 val , r 2 DFT ) used to calibrate the corresponding L n 1 and L n 2 norm deviations of the Fouriersynthesized distributions (S n,val , S n,tot and r 2 S n,val , r 2 S n,tot ). m = 1 m = 2 L m (r u.c. , val ) 0.04140798 0.00002884 L m (r u.c. , DFT ) 0.10471686 0.00204269 L m (r val , val ) 0.04533321 0.00003228 L m (r val , DFT ) 0.04700689 0.00003365 L m (r u.c. , r 2 val ) 0.10144559 0.00433510 L m (r u.c. , r 2 DFT ) 3.59309631 0.23552473 L m (r val , r 2 val ) 0.05329431 0.00004829 L m (r val , r 2 DFT ) 0.08058963 0.00014996 12 quantities are also calculated for the synthesized ED Laplacian distributions: in the formulas replace 'S n ' by 'r 2 S n ' and '' by 'r 2 '. Further criteria for the quality of the Fourier-synthesized EDs S n (r) are (i) the absence of unphysical negative ED values, and (ii) the absence of NNMs in agreement with the reference ED. The criterion of non-negative ED values was found to be violated by all non-weighted ('raw'), i.e. nonsmoothed, Fourier-synthesized all-EDs, a few non-weighted Fourier-synthesized val-EDs and a few Lan1D-synthesized all-EDs with a much too low weighting exponent p = 1.0. Application of smoothed Fourier synthesis corresponds in the first step to getting rid of negative EDs. But, as will be shown below, the requirement of the absence of artefact NNMs is a more severe challenge. It plays a major role in evaluating the synthesized ED and ED Laplacian distributions and adjusting the weighting-function exponent accordingly. Roughly speaking, for each Fourier-synthesis method presented herein, the appropriate exponent p of the weighting function p (H j ) was commonly adjusted (i.e. not for each resolution separately) to yield distributions without NNM artefacts for most resolutions, predominantly the higher ones. The choice to keep the exponent constant for each method was made to allow for studying systematic changes with increasing resolution.

Results: ED and ED Laplacian distributions from the method(exponent) type of Fourier synthesis
A look at the structure-factor amplitudes of the all-ED Fourier expansion reveals a very slow decay with increasing resolution (Fig. 3). From this figure, one would not have a criterion for determining a resolution necessary to yield a sufficiently converged Fourier-synthesized all-electron distribution. The valence-electron structure factors display a local minimum at about [sin()/] max ' 0.75 Å À1 , a local maximum at about [sin()/] max ' 1.3 Å À1 and an overall decay that is faster than the all-ED one. Note that the difference between the all-electron and the valence-electron structure factors is the core-electron structure factors (not shown). For the val-ED structure factors, the decrease of one order of magnitude from the local maximum at $0.75 Å À1 is only achieved at research papers Acta Cryst. Decay of valence-electron (blue) and all-electron (red) structure-factor amplitudes |F(H)| for CaB 6 with increasing resolution.

Figure 4
Decay of types of weighting function p used, shown for [sin()/] max = 1.5 Å À1 : (a) principal decay, p = 1; (b) smoothing used with various exponents 1.0 p 3.0; (c) complete prefactor (2H j ) 2 p (H j ) used for ED Laplacian synthesis. To distinguish the Lan1D smoothing curves for all-ED synthesis (2.25 p 3.0) from the val-ED ones (p = 1.0), they are always shown as pure lines skipping the data points. $3.25 Å À1 , which would be too high for experimental studies. The situation is even worse for the all-ED structure factors.
For raw Fourier synthesis, the similar size of neighbouring structure factors increases the possibility of sudden topological changes on inclusion of further neighbours, i.e. a nonconvergence of topological features, which is the typically observed behaviour. This can be prevented using smoothing factors p (H j ) as described above. In the following, the raw Fourier synthesis, which corresponds to constant weights p (H j ) = 1 for all structure factors F(H j ), will be considered the most aggressive weighting scheme compared with all other ones with smaller and variable weights. As will be seen below, raw Fourier synthesis with constant weights [equation (3)] has a strong tendency towards creation of artificial NNMs, and even regions with negative ED values upon Fourier-series truncation.
With exponent p = 1, overall concave decay (negative curvature) of p (H j ) is found for all 1D Fejer-type methods, except Fej_stl with linear decay throughout the whole range [0, [sin()/] max ] [ Fig. 4(a)]. In contrast, the 3D Fejer-type Pep3D factors on average display an overall convex decay. The Lanczos-derived schemes yield an initially concave decay, but exhibit an inflection point towards the end of the interval and become slightly convex. The most aggressive weighting schemes are found to be the Fejer-type ones Fej_pcl and Fej_cnt, which display the highest weights for most of the resolution interval [0, [sin()/] max ]. In particular, for these methods, NNMs were found for several resolutions (Table 3). This behaviour was counteracted by increasing the exponent p until a reasonable amount of resolutions became free of NNM artefacts at 1 2 H max ! 1.1 Å À1 . This was the case at values p = 1.5 and 1.75 for Fej_cnt and Fej_pcl, respectively [ Fig. 4(b)]. The less aggressive Lanczos-derived weights did not yield NNMs for 1 2 H max ! 1.1 Å À1 at p = 1.0 for the val-ED structure factors, but for the smoothing of all-electron structure-factor synthesized all-EDs with the Lan1D method, p = 2.25, 2.50, 2.75 and 3.0 had to be employed. In general, exponentiation with p > 1.0 increases the initial decay of p (H j ) and introduces a concave tail region or increases the degree of concavity of the tail region.
With the initially ( 1 2 H max ' 0) convex and tailing concave shape, the sigma factors p (H j ) are similar to the Debye-Waller factor which decreases according to exp(À 1 4 B iso H j 2 ). It is also interesting to see how the completely different methods Lan1D and Fej_cnt showing clearly different weighting factors p for p = 1.0 [ Fig. 4(a)] accidentally obtain very similar weights upon using the finally adjusted Fej_cnt exponent of p = 1.5 [ Fig. 4(b)]. From a signal processing point of view, the p factors for Fourier synthesis of the ED act as a low-pass filter. For the Fourier synthesis of smoothed ED Laplacian distributions, the p factors play a decisive role by over-compensating the quadratically increasing prefactor (2H j ) 2 , which would lead to non-convergence of the raw Fourier series. For the ED Laplacian reconstruction, the p factors act as a bandpass filter [ Fig. 4(c)].

Fourier synthesis of val-ED and val-ED Laplacian distributions
Fourier synthesis [equations (15), (16)] of val-ED S n,val (r) and val-ED Laplacian distributions r 2 S n,val (r) employs static val-ED structure factors. In view of the intended QTAIM analysis, the derived [equations (19), (20)] all-electron distributions S n;val ðrÞ and r 2 S n;val ðrÞ with added core distributions of the reference DFT wavefunction are always analysed in the following. The difference between these all-electron distributions and the corresponding reference DFT distribution is equal to the difference between the corresponding val-ED and val-ED Laplacian distributions [equations (21), (22)].
Raw Fourier synthesis of val-ED without smoothing [equation (3)] results in NNMs even at resolution 1 2 H max = 1.5 Å À1 [ Fig. 2(b)]. Therefore, for each method employed the exponent p of the smoothing factor p (H) was increased stepwise until a suitable number of resolutions in the range 0.5 Å À1 1 2 H max 5.0 Å À1 were free of NNMs. In order to study the important systematic changes for each method, the exponent p was not adjusted for each resolution separately. As a result, for each method a number of NNMs occurring at lower resolutions were accepted (Table 3).
Analysis of val-ED norm deviations. Raw Fourier synthesis of val-ED distributions is expected to converge in the L n 2 norm, which is depicted in the left panel of Fig. 5(b) Table 3 Fourier synthesis of val-ED S n,val (r) for CaB 6 .
Location of NNMs (see column 'location') and QTAIM basin populations from S n;val ðrÞ for all methods (exponents p) and resolutions 1 2 Hmax used. Method ones for the smoothed Fourier syntheses. The opposite behaviour is found for the L n 1 deviations at higher resolutions [except Pep3D(2.0), Fej_stl(1.0)], where the raw synthesis even displays a local maximum at about 1.2 Å À1 [ Fig. 5(a), left panel]. The different behaviour of the raw and smoothed ED distributions indicated by L n 1 and L n 2 norms is a result of the raw Fourier synthesis improving more on the higher deviations than the smoothed ones. This can be seen from the maximum norms given in the left panel of Fig. 5(c). The monotonic decay of the L n 1 norms for the smoothed Fourier syntheses is the mathematically expected behaviour. It is related to the expected uniform convergence for these weighting functions methods(p) as shown by the monotonic decay of the corresponding maximum norms in the left panel of Fig. 5(c). The overall decay of the maximum norm for the raw Fourier synthesis and the convergence indicated at higher resolutions indicate an unexpected uniform convergence of the raw val-ED Fourier synthesis in the full region.
The norm deviations in the valence region r val , shown in the right panels of Figs. 5(a)-5(c), are always smaller than the corresponding ones in the total region r u.c. [Figs. 5(a)-5(c), left panels], but to a variable degree. Notably, the uniform convergence for the smoothed synthesis is also found to be valid for the valence region [Figs. 5(a), 5(c), right panels], which was not initially expected.  In the valence region, the deviations L n 1 , L n 2 and Á max [Figs. 6(a)-6(c), right panels] of the synthesized val-ED Laplacian feature partially non-monotonic decays at lower resolutions and uniform convergence at higher resolutions for all functions [very slowly for Pep3D(2.0)]. The observed monotonic decays of the L n 1 norms in the valence region found for all smoothing functions besides the Fej_pcl(1.75) one may also suggest a certain regular convergence for the critical point Laplacian values. However, as can be seen in the critical points section below, all methods (with their specific exponents p used here) besides Pep3D(2.0) and Lan3D(1.0) display a similar damped oscillatory convergence behaviour with respect to deviation from the reference critical point values.

Analysis of val-ED
QTAIM basin analysis of S n;val ðrÞ. Classical QTAIM topological analysis was performed on the all-ED S n;val ðrÞ obtained by adding the correct core-ED of the reference DFT calculation to the synthesized val-ED [equation (19)  of the total electron number 50 e À /u.c. showed an error of maximally 0.0002 electrons in all investigations. The occurrence of NNMs for some resolutions (Table 2) does not lead to visible effects in the basin population and volume curves shown (Fig. 7).
All weighting functions methods(p) investigated display a rather regular convergence of the QTAIM atomic charges with increasing resolution [ Fig. 7(a)]. For all functions except Fej_stl(1.0) and Pep3D(2.0), a chemically reasonable accuracy was obtained at resolutions beyond 0.75 Å À1 . With respect to the decay of the weighting factors [ Fig. 4(b)], these functions correspond to the least aggressive ones. They are found to display the larger ED norm deviations (Fig. 5). Interestingly, the Ca atomic volumes [ Fig. 7(b)] do not simply display a similar behaviour to the charge curves, which is most clearly seen for function Fej_stl(1.0). Despite the comparably bad charge reconstruction, its volume reconstruction is comparable with the more competitive functions.
Critical point analysis of S n;val ðrÞ. The most challenging quality aspect for the Fourier-synthesized all-ED S n;val ðrÞ [equation (19)] and all-ED Laplacian r 2 S n;val ðrÞ [equation (20)] distributions is their ability to reproduce the chemical content contained in the original DFT-calculated distribution. For this purpose, critical points for all Fourier-synthesized ED distributions have been determined and analysed. The chemical focus requires (based on the DFT ED and chemical bonding arguments) the increase of ED at critical points in the sequence À 8 (c.c.p.-B 6 ) < À 6 (r.c.p.-BBB) < À 4 (b.c.p.-BB endo ) < À 5 (b.c.p.-BB exo ), negative values of the ED Laplacian for À 4 (b.c.p.-BB endo ), À 5 (b.c.p.-BB exo ), À 6 (r.c.p.-BBB), and positive ones for À 3 (b.c.p.-B 3 Ca) and À 8 (c.c.p.-B 6 ) increasing according to À 5 < À 4 < À 6 < À 3 < À 8 (cf. Fig. 1). These relations of signs and values are found to already be obeyed at lower resolutions for all Fourier-synthesis weighting functions methods(p) except for Pep3D(2.0) (Figs. 8, 9), even in those cases where at the original position of À 5 (b.c.p.-BB exo ) an NNM has been found. A rather monotonic approach of the corresponding ED values towards the initial ones is observed for all these functions (Figs. 8,9,left  Another important point concerns the size of the deviations at each critical point r c.p. for a given resolution. It is already obvious from the figures (Figs. 8, 9) that the absolute deviations of the synthesized ED and ED Laplacian values with respect to the reference values at each type of critical point are rather different. This observation is also valid for the relative deviations, given by Ár 2 S n;val ðr c:p: Þ rel: ¼ r 2 S n;val ðr c:p: Þ À r 2 DFT ðr c:p: Þ Â Ã = r 2 DFT ðr c:p: Þ : As an example, a collection of absolute and relative deviations from Lan1D(1.0, 2.0) Fourier synthesis are given for resolutions of 1.1 and 1.3 Å À1 in Table 4. The results for Lan1D(1.0) at 1 2 H max = 1.1 Å À1 could be considered as the lowresolution limit of a successful Fourier synthesis with method Lan1D being free from NNMs. However, within this series of Fourier syntheses the oscillatory approach to the reference values and visual ED Laplacian contour ripples (see the next section) for the higher resolutions (not yet for this comparably low resolution) may be disturbing. Certainly, the exceptionally high relative deviation of the ED Laplacian at r.c.p.-BBB of +87.7% is connected with this feature. A deviation of +100% would make the Laplacian value equal to 0, which would be a failure to qualitatively reproduce the chemical bonding condition r 2 (r.c.p.-BBB) < 0. An alternative candidate is Lan1D(2.0) at 1 2 H max = 1.3 Å À1 , where a lower resolution cannot be chosen because of NNM occurrence. It has the advantage of being a member of the Lan1D(2.0) series with a more regular convergence behaviour and negligible ED   (Table 3).
Laplacian ripples (see the next section) not only for this resolution, but also for higher ones up to 5.0 Å À1 . Obviously, there is not one best choice of low-resolution limit, it rather depends on priorities set up by the requirements of the investigation. This is also valid for the other methods presented (Figs. 8, 9).
All-ED S n;val ðrÞ and all-ED Laplacian r 2 S n;val ðrÞ distributions. In order to get a visual impression for the all-ED S n;val ðrÞ and all-ED Laplacian distributions r 2 S n;val ðrÞ obtained [equations (19), (20)] from Fourier-synthesized val-ED S n,val (r) and val-ED Laplacian r 2 S n,val (r), a selection of corresponding all-ED and all-ED Laplacian maps is presented in Figs. 10 and 11, respectively. For comparison with the Fourier-synthesized all-ED and all-ED Laplacian distributions (see below), those obtained with method Lan1D have been chosen here as well. For this method the value p = 1.0 of the weighting-function exponent could be selected for val-ED synthesis, which led to non-occurrence of NNMs for 1 2 H max ! 1.1 Å À1 and for 1 2 H max = 0.5 Å À1 . All synthesized ED distributions are completely smooth (Fig. 10), in contrast to the more challenging ED Laplacian distributions (Fig. 11 (Table 3). column). It can be seen that increasing resolution 1 2 H max is accompanied by increasing amount of contour ripples in the ED Laplacian distributions, which is a consequence of the rather aggressive weighting (i.e. comparably large weights).
For illustration of the effect of less aggressive weighting (smaller weights) on the Laplacian distribution smoothness, a comparison with a higher smoothing exponent p = 2.75 is shown in Fig. 11 Resolution-dependent all-ED Laplacian distributions r 2 S n;val ðrÞ in CaB 6 from val-ED Laplacian synthesis [equation (20)] for each resolution using method Lan1D. (a)-(c) employing exponent p = 1.0; (d)-(f) p = 2.75. The molecular graphs with b.c.p.'s (green), r.c.p.'s (blue), and c.c.p.'s (black) are shown; isolines are drawn at intervals of 0.05 e À bohr À5 . Note that the topology of (d) is substantially different from all the other ones: there is an NNM at the octahedron centre, and b.c.p.-BB endo and r.c.p.-B 3 are missing.

Figure 10
Resolution-dependent all-ED distributions S n;val ðrÞ in CaB 6 from val-ED synthesis [equation (19)] for each resolution indicated using method Lan1D. (a)-(c) employing exponent p = 1.0; (d)-(f) p = 2.75. The molecular graphs with b.c.p.'s (green), r.c.p.'s (blue) and c.c.p.'s (black) are shown; isolines are drawn at intervals of 0.01 e À bohr À3 . Note that the topology of (d) is substantially different from all the other ones: there is an NNM at the octahedron centre, and b.c.p.(B-B endo ) and r.c.p.(B 3 ) are missing. Table 4 Deviations of Lan1D(p) synthesized val-EDs and all-EDs and their Laplacians from reference DFT values at chemically important critical points for selected resolutions and exponents p.
Values given in the first row correspond to the reference DFT values (r c.p. ) and r 2 ðr c:p: Þ at each critical point; values given in each cell of the table body correspond to absolute and relative deviations [equations (44), (45)] of the ED (top row) and the ED Laplacian (bottom row), respectively; ED values and differences given in e À bohr À3 , ED Laplacian ones in e À bohr À5 . DFT  the ED Laplacian contour ripples can be smoothed out. No systematic efforts have been undertaken to find for each resolution the smallest necessary exponent p value for obtaining a certain predefined smoothness of the ED Laplacian contours. It may be sufficient to mention that already with p = 2.0 all ED Laplacian contours calculated with this model are visually smooth (a numerical criterion has not been employed). The value p = 2.75 was chosen here to make a visual comparison with the Fourier-synthesized all-ED distributions shown below, which have been obtained with the same exponent. For the synthesized val-ED distributions, the higher exponents p = 2.0 and 2.75 yield systematically higher overall L n 1 and L n 2 norm deviations of the val-ED in valence regions than with more aggressive weighting p = 1.0 [ Fig.  12(a)].
For the synthesized val-ED Laplacian the situation is more complex. The L n 1 (r val , r 2 S n,val ) deviations in the valence region [ Fig. 12(b)] are found to be smaller with aggressive weighting p = 1.0 than with p = 2.75 weighting until 1 2 H max = 2.0 Å À1 , and become larger for resolutions 1 2 H max ! 2.5 Å À1 .
For all resolutions shown [ Fig. 12(b)], the intermediate exponent p = 2.0 yields lower norm deviations L n 1,2 than the p = 2.75 one. Thus, for higher resolutions 1 2 H max ! 1.5 Å À1 , an exponent 1.0 < p 2.0 Å À1 might be sufficient to prevent visible Laplacian contour ripples. This may be related to the original Lanczos result (Lanczos, 1956) of p = 1 for gradient synthesis, implying p = 2 for ED Laplacian synthesis. It is, however, not clear how his results on purely 1D Fourier synthesis apply to a 1D condensed variant (like Lan1D) of a real 3D Fourier-synthesis task. In summary, it seems that the increasing rippling in the ED Laplacian distribution contours has an adverse effect on the L n 1 norms of the p = 1.0 syntheses, and that they improve upon smoothing out those ripples. Once the ripples are smoothed out (with p = 2.0), the norm deviations become worse again upon increased smoothing (p = 2.75).
Certainly, the smoothing of ED Laplacian contours via increased p values simultaneously increases the ED L n 1 deviations [seen also locally in critical point analysis, Fig. 9(a), left panel] of the already smooth ED distributions as just indicated, which means that a suitable compromise p value must be found, depending on the focus of the application.

Fourier synthesis of all-ED and all-ED Laplacian distributions
As was shown in Fig. 3, the all-ED structure factors decay very slowly with increasing resolution 1 2 H max , which certainly makes Fourier synthesis in the atomic core regions a very challenging task. However, the question arises: to what extent these effects adversely affect ED analysis in the valence region? Because of its quality in reproduction of the chemical bonding features with exponent p = 1.0 for val-ED Fourier synthesis, method Lan1D was chosen to investigate this issue. By systematic p-value variation it was found that for all-ED Fourier synthesis, exponents in the range 2.25 p 3.0 are useful to prevent negative ED and NNM artefacts for most resolutions used (Table 5). The preferred model for systematic evaluations was the Lan1D(2.75) one, which showed no such artefacts for all resolutions. Other smoothing variants in the indicated range were also investigated, because for certain resolutions, i.e. 1 2 H max = 1.5 Å À1 , their Fourier syntheses were CaB 6 : valence region relative norm deviations L n 1 (full lines) and L n 2 (dashed lines) for Lan1D synthesis of (a) val-ED S n,val (r) and all-ED S n,tot (r), and (b) val-ED Laplacian r 2 S n,val (r) and all-ED Laplacian r 2 S n,tot (r). Val-ED and val-ED Laplacian synthesis with exponents p = 1.0 (black lines), 2.0 (blue lines), 2.75 (red lines), and all-ED and all-ED Laplacian synthesis with p = 2.75 (green lines).

Table 5
Fourier synthesis of all-ED. Analysis of all-ED norm deviations. Raw Fourier synthesis of the all-ED is expected to converge in the L n 2 norm, which is verified in the left panel of Fig. 13(b). The L n 2 deviations are smaller than the ones for the smoothed Fourier syntheses, and the convergence is faster. The opposite behaviour is found for the L n 1 deviations depicted in the left panel of Fig. 13(a). Here, the smoothed all-ED distributions have lower deviations than the raw one, which is a result of the raw Fourier synthesis improving more on the higher deviations occurring at higher values [ Fig. 13(c), left panel], namely in the core regions. The raw Fourier-synthesis L n 1 deviations display the occurrence of a local maximum at about 0.9 Å À1 and convergence at higher resolutions. The convergence of the L n 1 norm for the smoothed Fourier syntheses [ Fig. 13(a), left panel] is the mathematically expected behaviour. It is related to the uniform convergence obtained using these weighting functions methods(p) as can be seen by the monotonic decay of the corresponding maximum norms in the left panel of Fig.  13(c). The monotonic decay of the maximum norm on the grid for the raw Fourier synthesis indicates uniform convergence as well, which is surprising and the detailed reason is not known. Certainly, for all resolutions investigated, the maximum norm is found at the position of the Ca nucleus, which is located on a grid point. This can be considered a lucky circumstance, because it is the position where the maximum norm over the whole unit-cell region is expected to be located, mainly because the boron nucleus (not being located on the evalua- Fourier synthesis of all-ED for CaB 6 using method Lan1D with exponents p = 2.25, 2.50, 2.75 and 3.0. (a)-(c) Norm deviations L n 1 , L n 2 and Á max , respectively, of all-ED S n,tot (r); left panel, total u.c. region r u.c. ; right panel, valence region r val . tion grid) has a lower atomic number. It is important to realize that with large deviations from the reference DFT ED at the Ca nuclear position and their slow convergence indicated by values of DFT (r Ca ) ' 7522 bohr À3 and S n,tot (r Ca ) ' 713 bohr À3 for raw synthesis and 148 bohr À3 for Lan1D(2.75) at 1 2 H max = 5.0 Å À1 , Fourier synthesis of the all-ED in the core region near the nucleus is not suitable in a quantitative way with the method(exponent) type of approach applied here. The question remains as to whether all-ED Fourier synthesis in the valence region can give chemically meaningful results.
In the valence region, all deviation types L n 1 , L n 2 and Á max [Figs. 13(a)-13(c), right panels] for the smoothed all-ED distributions are much smaller than for the corresponding raw ones. While the maximum norm of the raw distribution displays an irregular behaviour up to high resolutions, the smoothed ones decay monotonically beyond a local maximum at about 1.2 Å À1 . It is noteworthy that the locations of the maximum norm in the valence region are always found at grid points directly in the neighbourhood of the exclusion radius of the boron atoms, such that their values are considered to be (and actually found to be, see below) significantly larger than the deviations expected at the critical point locations deep inside the valence region. Although overall monotonic decay of L n 1 and L n 2 in the valence region is not seen for the smoothed distributions, uniform convergence of Fourier-synthesized all-EDs is indicated by the norms L n 1 , L n 2 and Á max at higher resolutions [Figs. 13(a)-13(c), right panels]. The more aggressive weighting functions (lower p values, higher weights) yield the smaller deviations at these resolutions. Fourier synthesis of all-ED Laplacian for CaB 6 using method Lan1D with exponents p = 2.25, 2.50, 2.75 and 3.0. (a)-(c) Norm deviations L n 1 , L n 2 and Á max , respectively, of all-ED Laplacian r 2 S n,tot (r); left panel, total u.c. region r u.c. ; right panel, valence region r val .
Analysis of all-ED Laplacian norm deviations. For raw Fourier synthesis the all-ED Laplacian deviations in the full unit-cell region are divergent [Figs. 14(a)-14(c), left panel] as mathematically expected. While this is clearly indicated in the L n 1 and L n 2 norm deviations [Figs. 14(a), 14(b), left panel], convergent behaviour of the maximum norm Á max is simulated at lower resolutions until a sudden strong increase beyond 3.0 Å À1 changes the convergence scenario [ Fig. 14(c), left panel]. The abrupt rise is an artefact caused by the combination of two factors, namely (i) the finite grid resolution, and (ii) the electron-nuclear cusp (Kato, 1957). Due to the shape of the plane waves, the electron-nuclear cusp cannot be represented by Fourier synthesis with structure factors available at any experimental conditions and not in the present theoretical study either. In contrast, the electronnuclear cusp is present in the given reference DFT all-ED of CaB 6 and leads to a value of 0 for the DFT all-ED Laplacian given at the nuclear positions (in fact, the value is mathematically undefined, and is arbitrarily given a value of 0), such that these positions were excluded from the maximum norm evaluations of the synthesized all-ED Laplacian distributions. On the regular grid used, this only concerns the Ca position at (0, 0, 0), which is a grid point, while the B-atom positions are not located on grid points. Thus, on evaluation of the full unitcell maximum norm of the all-ED Laplacian, the highest deviation is found for a grid point closest to the B nucleus for all resolutions 1 2 H max 3.0 Å À1 . For all higher resolutions investigated, another grid point features the highest deviation. It is located one grid point away from the Ca position, i.e. at a larger distance from the Ca nucleus than the previous grid point was from the B nucleus. Principally, the region around the Ca nucleus is the one where all the highest deviations are expected to be located. However, due to the larger distance of the Ca nucleus to the next grid point, the point with the shorter distance to the B nucleus displays the highest deviation at lower resolutions, until the divergence of the large deviations close to the Ca nucleus finally dominates the maximum norms at higher resolutions. This nicely illustrates the difficulties that can arise on using the maximum norm as a sole criterion of convergence behaviour.
For the smoothed distributions in the total region a uniform convergence seems to be indicated by the    QTAIM basin analysis of S n;tot ðrÞ. Standard QTAIM topological analysis was performed on the all-ED distribution obtained by Fourier synthesis of the all-ED. The numerical electron recovery of the total electron number 50 e À /u.c. in all investigations showed an error of maximally 0.0002 electrons. The occurrence of NNMs for some resolutions (Table 2) does not lead to visible effects in the effective charge and volume curves shown (Fig. 15). The Ca effective charge obtained from QTAIM analysis of the reference DFT ED amounts to Q eff (Ca) = +1.52.
In contrast to the rather monotonic improvement of the Ca effective charges obtained from val-ED Fourier synthesis [ Fig.  15(a), Lan1D(1.0) val-ED results shown by the black line for comparison], all-ED synthesis features comparably large deviations improving in an oscillatory manner to the correct value [ Fig. 15(a)]. The more aggressive variants of Lan1D (lower p values) Fourier synthesis are found to display the larger oscillations. They show smaller L 1 norm deviations [ Fig. 13(a)] than the less aggressive ones with higher p values. An accuracy with respect to the Ca effective charge comparable with that found already at resolutions of about 0.75 Å À1 for the val-ED Fourier synthesis can be obtained with the all-ED synthesis only at resolutions beyond 2.5 Å À1 [ Fig. 15(a)]. While the Ca volumes are found to be systematically too small at all resolutions for the val-ED synthesis [ Fig. 7(b), Fig. 15 Resolution-dependent all-ED distributions S n;tot ðrÞ from all-ED synthesis for CaB 6 using weighting function Lan1D(2.75). The molecular graphs with b.c.p.'s (green), r.c.p.'s (blue) and c.c.p.'s (black) are shown for each resolution; isolines are drawn at intervals of 0.01 e À bohr À3 .

Figure 18
Resolution-dependent all-ED Laplacian distributions r 2 S n;tot ðrÞ from all-ED Laplacian synthesis for CaB 6 using weighting function Lan1D(2.75). The molecular graphs with b.c.p.'s (green), r.c.p.'s (blue) and c.c.p.'s (black) are shown for each resolution; isolines are drawn at intervals of 0.05 e À bohr À5 .
all-ED synthesis distributions as well, À 8 (c.c.p.-B 6 ) < À 6 (r.c.p.-BBB) < À 4 (b.c.p.-BB endo ) < À 5 (b.c.p.-BB exo ). The convergence of the ED deviations from the original values is non-monotonic for the three critical points À 5 , À 4 and À 6 with the largest ED values, displaying one local minimum at roughly 1 2 H max 1.1 Å À1 , but from this point on the deviations are found to decrease monotonically with increasing resolution.
The absolute and relative deviations of the synthesized all-ED and ED Laplacian values are given in Table 4 for resolutions 1 2 H max = 1.3 and 1.5 Å À1 . Inspecting the deviations for all-ED Laplacian synthesis at 1 2 H max = 1.3 Å À1 , the relative deviation of +144% at r.c.p.-BBB reveals that the chemical bonding criterion of a negative all-ED Laplacian is even missed there. This is no longer the case for resolutions !1.5 Å À1 . While the absolute values of the ED deviations for Lan1D(2.75) all-ED synthesis at 1.5 Å À1 are of roughly similar size as those for Lan1D(2.0) val-ED synthesis at 1 2 H max = 1.3 Å À1 , the corresponding all-ED Laplacian deviations for the Lan1D(2.75) all-ED synthesis are significantly higher at the three critical points, with negative reference ED Laplacian values compared with Lan1D(2.0) val-ED Laplacian synthesis. Similar to the situation discussed above for val-ED Laplacian synthesis with Lan1D(1.0), also for all-ED Laplacian synthesis using Lan1D(2.75) an oscillatory approach to the reference value is found [ Fig. 16(b)] for critical points in the resolution range 1.0 Å À1 1 2 H max 2.0 Å À1 , while for resolutions beyond 2.0 Å À1 a rather monotonic decrease of the deviations from the reference DFT values can be expected. Also, for all-ED Laplacian synthesis, increasing the p value of Lan1D(2.75) is expected to smooth out this oscillatory approach at the price of increased deviations for the synthesized all-ED values. This has not been explicitly tried here. For the present study, it is considered to be sufficient to have established a rough comparison between valence-and all-ED Fourier-synthesis efficiency (Table 4) for the example of CaB 6 and choice of synthesis method Lan1D. A wider comparison will be given in Section 4.
All-ED S n;tot ðrÞ and all-ED Laplacian r 2 S n;tot ðrÞ distributions. For visual inspection, the synthesized all-ED and all-ED Laplacian distributions are depicted in Figs. 17 and 18, respectively. It can be seen that not only the all-ED contours, but also the corresponding Laplacian contours, are very smooth. Thus, the comparably high weighting function exponent p = 2.75 dictated by avoidance of artefact NNMs seems to be sufficient for smoothing out Fourier-synthesis ripples as well. This is different from the val-ED synthesis situation, where the avoidance of NNMs was already achieved at p = 1.0, but val-ED Laplacian contour ripples at higher resolutions were clearly visible [ Figs. 11(b), 11(c)]. They were seen to be smoothed out at the same exponent p = 2.75 as well [Figs. 11(e), 11(f)]. Nevertheless, the norm deviations for all-ED and all-ED Laplacian synthesis in the valence regions are found to be significantly higher than for val-ED [ Fig. 12(a)] and val-ED Laplacian synthesis [ Fig. 12(b)] using the same exponent p = 2.75.

Discussion
Six different weighting function methods, namely Lan3D, Pep3D, Lan1D, Fej_pcl, Fej_cnt and Fej_stl, for weighted Fourier synthesis [equations (15), (16)] have been employed for investigation of the reconstruction quality of characteristic all-ED and all-ED Laplacian features found in the corresponding reference DFT distributions of CaB 6 . Two kinds of all-ED and all-ED Laplacian reconstruction procedures have been considered: (i) Fourier synthesis of val-ED S n,val (r) and val-ED Laplacian r 2 S n;val with subsequent addition [equations (19), (20)] of the corresponding reference DFT core distributions resulting in distributions S n;val ðrÞ and r 2 S n;val ðrÞ, respectively, and (ii) Fourier synthesis of all-ED S n,tot (r) and all-ED Laplacian r 2 S n;tot ðrÞ distributions.
The actual weighting functions employed are characterized by the combination of a method with a smoothing exponent p, denoted as method(p). It is noteworthy that, while the nonexponentiated method functions method(1) (i.e. p = 1) display rather different weighting curves [ Fig. 4(a)], the combination of each method with a different exponent can make them more similar. For example, weighting function values of Fej_cnt(1.5) turned out to be very similar to the ones for Lan1D(1.0) [ Fig. 4(b)], which was initially not the case [ Fig.  4(a)]. Nevertheless, although overall quality, as measured by norm deviations, and chemical quality, as measured from specific QTAIM-derived deviations, were found to be very similar for both functions above, the tiny differences  Table 6 Summary of numerical quality criteria fulfilment based on QTAIM analysis for various weighting functions method(exponent p) tested.
Values are given for val-ED synthesis ('S n,val synthesis') and all-ED synthesis ('S n,tot synthesis'). Resolutions for criterion C{r 2 (r c.p. )} are put into brackets, where the ED Laplacian at critical points does not clearly converge. For each method(exponent) the decisive resolution for fulfilment of all four criteria is marked in bold. remaining led to non-occurrence of NNM artefacts at slightly different resolutions (Table 6). Norm deviations of the raw and synthesized ED and ED Laplacian distributions in the whole unit-cell region have been computed to verify the mathematically expected convergence behaviours with increasing resolution. Furthermore, as an intermediate step towards the more point-wise convergence investigation at selected ED critical points in the QTAIM framework, norm deviations in the valence region of the unit cell have been studied. Considering the intended QTAIMtype topological analysis, the observed convergence of the L n 1 norm valence-electron and all-ED deviations in the valence region [Figs. 5(a) and 13(a), right columns] is an important reconstruction quality indicator, because interatomic separatrices (basin surfaces) and critical points located on them are then embedded in a region with systematic improvement of the overall linear deviations L n 1 (homogeneous convergence). The same is valid for the more sensible ED Laplacian distributions [Figs. 6(a) and 14(a), right columns].

Resolutions
Finally, the convergence of QTAIM analysis based results with increasing resolution was evaluated for each weighting function method(exponent p). The quantities considered are the effective atomic charges, and ED and ED Laplacian values at chemically important critical points. In the case of val-ED synthesis with weighting functions Lan1D(1.0), Fej_pcl(1.75), Fej_cnt(1.5) and Fej_stl(1.0) convergent behaviour is detected (Figs. 7-9). The val-ED Laplacian displays a kind of damped oscillatory approach to the reference DFT values, and the weighting functions Lan3D(1.0) and Pep3D(2.0) even showed some difficulties converging the val-ED Laplacian at selected critical points (Fig. 8) for larger resolutions (e.g. 5 Å À1 ), such that their occurrence in the final evaluation of reconstruction quality is only chosen for the sake of completeness, and the problematic non-converged criterion is put into brackets in Table 6. This problem of the Lan3D(1.0) scheme needs further investigation in the future.
For all-ED and all-ED Laplacian synthesis using weighting function Lan1D(2.75), convergence of the QTAIM effective charge Q eff (Ca) was found to show oscillatory behaviour with much larger deviations than found for val-ED Lan1D(1.0) synthesis at the same resolution (Fig. 15). The reason is the large average all-ED set up by structure factor F tot (000) of 50 e À /V u.c. , which is compared with the average val-ED of 20 e À /V u.c. from F val (000). This constant all-ED level must be redistributed by the subsequent structure factors to build up the comparably huge and narrow atomic 'peak' in the core regions and depleting the valence regions between the atoms accordingly. The redistribution process works rather slowly with increasing resolution, as has been indicated above from the synthesized values of the all-ED at the Ca position at the resolution of 5 Å À1 to be only 9.5% (raw synthesis) and 2% [Lan1D(2.75)] of the reference DFT value. As a consequence, the nuclear core 'peak' is already too wide from using raw Fourier synthesis [equations (15), (16), with constant p (H j ) = 1] and even wider for smoothed Fourier synthesis [equations (15), (16), within the method(exponent) type of approach], leading to decreased charge transfer for lower resolutions, where these core-electron regions significantly overlap. For all-ED synthesis with resolutions !2.0 Å À1 Q eff (Ca) is found to become reasonably converged (Table 6), indicating reasonably reduced core(Ca)-core(B) ED overlap. The synthesized all-ED values at the characteristic critical points investigated show rather smooth convergence behaviour for resolutions !1.1 Å À1 , while all-ED Laplacian values again (like for the val-ED Laplacian) show an oscillatory approach to the reference DFT values (Fig. 16).
In a specific comparison of the reconstruction qualities of each weighting function method(p) used, a collection of four criteria has been set up; classification of the most successful function is the one fulfilling all criteria at the lowest resolution. The criteria chosen are based on chemical bonding arguments in the framework of QTAIM, and also take into account the precision necessary for typical chemical bonding discussions. They are clearly not universally objective, but serve for the predefined purpose.
The first and most important criterion for comparison of the reconstruction quality is the absence of artefact NNMs in the synthesized EDs. As already mentioned, this criterion was used to select the lowest possible smoothing exponent p ! 1.0 for each weighting method.
A second criterion is based on the reconstruction of the effective charges Q eff of the different species, i.e. reconstruction of Q eff (Ca) = +1.52 obtained from QTAIM analysis of reference DFT ED is sufficient. A satisfactory recovery of the reference Ca-atom atomic charge is considered to be achieved within a deviation |ÁQ eff (Ca)| 0.05, which is fulfilled for most val-ED Fourier-synthesis weighting functions of type method(exponent) at resolutions ! 0.75 Å À1 .
Based on the reference DFT ED and chemical bonding arguments, an increase of ED at critical points in the sequence À 8 (c.c.p.-B 6 ) < À 6 (r.c.p.-BBB) < À 4 (b.c.p.-BB endo ) < À 5 (b.c.p.-BB exo ) was required to be obeyed for a successful Fourier synthesis. This third evaluation criterion, denoted C{(r c.p. )}, is found to be obeyed at virtually all resolutions investigated (Table 6).
Satisfactory reconstruction of the val-ED Laplacian distribution is more challenging. One problem is the kind of damped oscillatory approach of the synthesized values towards the reference ones (Figs. 8,9) found for all weighting functions [besides Pep3D(2.0)]. As investigated for method Lan1D, a suitable increase of weighting function exponent p damps the oscillatory approach and can thus decrease the lowest resolution for fulfilment of condition C{r 2 (r c.p. )} for the other methods with comparably small exponents 1.0 p < 2.0. The complete fulfilment of the requests[requirement for] of negative values of the ED Laplacian for À 4 , À 5 , À 6 , and positive ones for À 8 and À 3 (b.c.p.-B 3 Ca) increasing according to À 5 < À 4 < À 6 < À 3 < À 8 (which is the fourth evaluation criterion denoted C{r 2 (r c.p. )}) is mainly dependent on the ED Laplacian sign at À 6 (r.c.p.-BBB). The reference value of À0.0279 bohr À5 is very close to zero, such that resolutiondependent, rather small variations in absolute value may yield positive ED Laplacian values here. Hence, the fulfilment of relations C{r 2 (r c.p. )} is strongly method and exponent dependent and is found to vary between resolutions of [sin()/] max ! 0.5 Å À1 [functions Fej_cnt(1.5) and Lan1D(1.0, 2.0)] and [sin()/] max ! 1.8 Å À1 [function Fej_pcl(1.75)].
The simultaneous fulfilment of all four conditions based on QTAIM analysis of the synthesized ED and ED Laplacian distributions of CaB 6 leads to the conclusion that for val-ED Fourier synthesis a resolution of at least 1.2 Å À1 is necessary for chemically suitable ED and ED Laplacian distributions. For comparable all-ED synthesis quality, a data set resolution of 2.0 Å À1 is necessary (Table 6). It is noteworthy that the high value for all-ED synthesis is caused by the |ÁQ eff (Ca)| condition and not the C{r 2 (r c.p. )} condition, which would have already been satisfied at and beyond 1.7 Å À1 .

Conclusions
The topic of this pilot study is the development and evaluation of a strategy for the task of extracting sufficiently precise ED and ED Laplacian distributions using a Fourier 'back-transformation' process with an incomplete number of structure factors (Fourier synthesis). The compound CaB 6 , crystallizing in a cubic cP7-type of structure, has been chosen for this study. As a result, a Fourier-synthesis approach has been presented, which yields a systematic reconstruction of ED and ED Laplacian distributions from quantum-chemically (DFT) calculated valence-electron and all-ED static structure factors of variable resolution. The decisive issue is the application of suitable weighting functions to avoid series termination artefacts, while extracting the maximum amount (precision) of chemical bonding information possible.
The features of the ED and ED Laplacian reconstructions obtained by the present Fourier-synthesis approach have been studied for six weighting methods, namely Pep3D, Lan3D, Lan1D, Fej_pcl, Fej_cnt and Fej_stl. Pep3D and Lan3D have been explicitly adopted from the literature, methods Fej_pcl, Fej_cnt and Fej_stl are newly developed, and method Lan1D has been mentioned before in the literature, but not explicitly formulated. For the purpose of getting rid of NNM artefacts, a new strategy has been introduced to supplement each weighting method with a suitable smoothing factor exponent, such that the final weighting functions used were of the type method(exponent) (ME approach). The exponent p was adjusted to take the smallest possible value (p ! 1.0) leading to ED distributions without NNM artefacts. The smaller p values led to smaller norm deviations of the synthesized EDs. Using higher exponents within each method, eventual resolution-dependent ED Laplacian contour ripples can be systematically smoothed out as well, although at the price of increased ED norm deviations, such that a compromise needs to be found.
Convergence of the ED and ED Laplacian distributions with respect to the corresponding reference DFT-based distributions with increasing resolution [sin()/] max was clearly demonstrated by analysis of norm deviations of the synthesized distributions, QTAIM effective charge deviations, and ED and ED Laplacian value deviations from the reference DFT values at chemically important critical points. The criteria for successful sufficiently precise reproduction of the characteristic ED and ED Laplacian features important for chemical bonding arguments were chosen to be rather low, e.g. they were based on qualitative relations of values between different critical points and not on absolute value convergence. Based on these types of criteria, successful reconstruction of the val-ED and its Laplacian from valenceelectron structure factors is found for resolutions !1.2 Å À1 and for all-ED structure factors !2.0 Å À1 . In the evaluations of the reconstruction quality, the 1D weighting methods have overall been more successful than the 3D ones, with the best results obtained from Lan1D, Fej_pcl, Fej_cnt. This list of suitable methods is not considered to be a static one; it may be extended by further methods in the future, and some of the present methods may turn out to be less suited for, e.g., systems with lower symmetry. Generally, the methods presented here are not suitable for a realistic reconstruction of the all-ED in the atomic core regions, because reconstruction of the steep nuclear all-ED 'peak' with increasing resolution is too slow. An approach complementary to the present one, focusing mainly on reproduction of the all-ED nuclear 'peaks', has been reported by Altomare et al. (2008).
The prospects of application of the presented ME type of Fourier-synthesis approach could be rather wide, provided some substantial supplementary studies are undertaken in the future.
In combination with experimental ED reconstruction studies based on the Hansen-Coppens model, the presented approach could be applied without further modification. The HC model would then play the role of the reference DFT model in the present study, and the fitted multipoles would correspond to the val-ED val (r), whose static structure factors are accessible. The core-ED in the HC model is derived from free atoms and takes the role of core (r) in the present study. Calculation of static ED distributions within the HC model using the multipole parameters corresponds to an extrapolation of the experimental data set to infinite resolution. Conversely, in the framework of Fourier synthesis, the application of mathematical weighting functions decaying with increasing sin()/ until [sin()/] max corresponds to systematically down-weighted contributions of the higher reflections, which could be considered as an under-interpretation of these data. Therefore, a less biased strategy for robust experimental static ED studies could be to consider both approaches in parallel. For gathering more experience with Fouriersynthesis techniques for ED and ED Laplacian reconstruction, further theoretical studies are necessary, e.g. for non-cubic structures.
The applicability of the presented method(exponent)-type approach to Fourier synthesis of dynamic ED and ED Laplacian distributions in the valence regions is expected. The corresponding studies could be done on the basis of the HC model using the core-electron and valence-electron structure factors and the thermal parameters. Because of the physical smoothing of the distributions caused by thermal averaging,